Density of States and Critical Behavior of the Coulomb Glass 
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We present zero-temperature simulations for the single-particle density of states of the Coulomb glass. Our 
results in three dimensions are consistent with the Efros and Shklovskii prediction for the density of states. 
Finite-temperature Monte Carlo simulations show no sign of a thermodynamic glass transition down to low 
temperatures, in disagreement with mean-field theory. Furthermore, the random-displacement formulation of 
the model undergoes a transition into a distorted Wigner crystal for a surprisingly broad range of the disorder 
strength. 
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The Coulomb glass (CG) is the earliest paradigm for under- 
standing the effects of strong disorder in electronic systems 
with long-ranged interactions. Among its applications are the 
space-energy correlations in transistors, the magnetization- 
switched metal-insulator transitions in tunnel devices, the 
cotunneling magnetoresistance in ferromagnetic devices, the 
ambipolar gate effect, the huge magnetoresistance in semi- 
conductor stacks, and the transparent refractory oxides, to 
name a few. In the CG, the Coulomb interaction remains 
long ranged because the disorder localizes the electrons and 
thus impedes screening. Therefore, the system forms its 
low energy states by long-range configurational changes and 
avalanches. After some early approaches Efros and 

Shklovskii (ES) argued that the stability of the low-energy 
states against the long-ranged single-particle dynamics re- 
quires the formation of a soft "Coulomb gap" in the single- 
particle density of states (DOS) of the form p{E) ^ \E\^ 
with (5 = 13 — 1 as an upper bound and where the energy 
E is measured from the Fermi level Istl; D being the space 
dimension. This Coulomb gap leads to a typical variable 
range hopping form of the low-temperature conductivity, i.e., 
cr(T) = croexp[-(ro/r)-i/2]; ctq and Tq constant. 

There is considerable experimental evidence in support 
of thesepredictions from transport measurements of a{T) 
10. S 0,01, as weU_as from tunneling conductance measure- 

Bin. 



ments of p{E) 



However, subsequent theoreti- 



cal considerations arrived at an exponential form of the DOS 
by considering multi-electron "polaronic" processes; p{E) ~ 
exp [—{Eq/EY/'^] 112111311 . throwing the status of theoretical 
predictions for the DOS in question. 

A large number of nonequilibrium glassy phenomena have 
been observed in disordered electronic systems. These in- 
clude slow dynamics lfl4l [Tsi Til 171. aging, and memory 
effects [l9[l2oll . as well as changes in the noise spectrum 
1 2111 . However, the existence of a thermodynamic glass tran- 
sition cannot be directly surmised from glassy dynamics. In 
fact, no well-defined thermodynamic glass transition has been 
found in association with these phenomena to date in three- 
dimensional (3D) systems. 

These different theoretical predictions merged with the in- 



sight of Pastor and Dobrosavljevic, who described a disor- 
dered electron system with long-range interactions u sing th e 
replica-based theoretical framework of glass physics l22ll23ll . 
Their work offered a unified platform to analyze both the DOS 
and the glassy characteristics of these systems. This pro- 
gramme was subsequently expanded by the work of Miiller, 
loffe and Pankov who included replica symmetry breaking 



technology into their calculations 1124 |25L 12611 . All these 
studies concluded that — within a mean-field approach — a soft 
Coulomb gap exists in the single-particle DOS at T = 0. Fur- 
thermore, for T < Tc - W-^l"^, where W is a measure of 
the disorder, the system freezes into a "Coulomb glass" state. 
Note that the Coulomb glass is analogous to a spin glass in a 
(random) field 12711 which is known to not order. 

The Coulomb glass has attracted considerable attention nu- 
merically as well. The initial work by Davies, Lee and 
Rice reported the observation of a soft gap, but the data 
were not conclusive with respect to the detailed functional 
form of the DOS 12811 . Subsequent numerical studies rep- 
resented the disorder either by random site energies (CG) 
HlHIIllllElfsil or by random displacements (RD) be- 



tween the sites ISllMlSllMlMIM]. While the CG and RD 
models have differeiit symmetries (and thus possibly different 
universality classes 14111 ). it has nevertheless been argued that 
they both ade quat ely capture the key aspects of the real elec- 
tronic system 13811 . Different studies of the DOS in 3D have 
reported a DOS vanishin g at the Fermi level with p(E) ~ \E\^ 
with ,5 = 2.1 - 2.6 mE. 1|]. Even more surprising was the 
claim of a strongly disorder-dependent exponent 5 13811 . Fur- 
thermore, studies attempting to locate a transition to a glassy 
state were only successful in the RD model 

The state of the field can be summarized as follows: A 
soft gap in the DOS has been widely confirmed, but the pre- 
dicted ES exponent (S = D — 1 is consistent only with exper- 
imental data, not with numerical simulations. A true finite- 
temperature transition to a glassy state has numerical support 
in the RD model but lacks evidence in the CG model and in 
experiments. 

Our results show that in the 3D CG model 5 is close to 2 and 
weakly disorder dependent lE6ll and we find no signature of 



2 



a finite-temperature glass transition. Furthermore, in the RD 
model the low-temperature ordering is indicative of a distorted 
Wigner crystal. 

Model and numerical details. — The Coulomb glass (CG) 
Hamiltonian is given by 



1 

^CG = ^ - — -^)+Yl 



(1) 



where € {0, 1} is the electron number at site i, v the fill- 
ing factor, and ep' / nrij the Coulomb repulsion. The sites lie 
on a three-dimensional lattice of size N ^ L^, and the elec- 
tron number is coupled to Gaussian-distributed random site 
energies Si with zero mean and standard deviation W, i.e.. 
Vie,) = (27riy2)-i/2 cxp(-e|/2W^2)^ In the RD model, 

instead of random site energies, the disorder is represented 
by Gaussian-distributed random displacements of the lattice 
sites with standard deviation ^/3W. The DOS is given by 
the disorder average of p{E) = (l/N) - E,) with 

-^i — YLj^ii'i^j ~ / + £j the local single-particle 
energy yj. 

For the simulations we use particle-conserving dynamics 
and periodic boundary conditions. To cope with the long- 
range Coulomb interactions we perform a resummation tech- 
nique in which we sum all interactions over periodic im- 
ages and renormalize the energy scales such that the nearest- 
neighbor distance is a = 1. To compute the ground-state 
DOS (T = 0) we use extremal optimization 14211 . For the 
CG model we perform updates and study systems of up 
to = \A'^ sites in 3D for W = 0.2 and 0.4 and average over 
3000 disorder samples for L<12 and 1800 (800) samples for 
L = 14 for = 0.2 {W = 0.4). For the RD model we study 
N = 14"^ sites and average over 100 disorder samples (fluc- 
tuations are small). For the study at finite temperatures we 
use exchange Monte Carlo 1I43I Hi . Equilibration is tested 
by a logarithmic data binning. Once the last three bins agree 
within errors, the system is in thermal equilibrium. Simula- 
tion parameters can be found in Table U 

Results for the density of states. — Figure[T|(top and cen- 
ter panels) show the DOS at T = for the 3D CG model 
for two disorder strengths close to the Fermi level {E = 0) 
at half filling {v = 1/2); the insets show the whole func- 
tional shape. The data can be fit very well with a form ~ |£'|'^ 
with 5 = 2.01(2) (L = 14) for W = 0.2 and 5 = 1.83(3) 
[L = 14) for W = 0.4 (restricted to \E\ < 0.3), which is 
close to the ES value of (5 « D — 1. 

Fig. [T] (bottom) shows the DOS of the RD model for L = 
14. The DOS shows a pronounced double-peak, the width of 
the peaks dependent on W. There is no sign of the character- 
istic Coulomb gap shape, moreover the peaks at \E\ ^ 1 are 
typical of a Wigner crystal (WC). Thus the DOS of the RD 
model is indicative of the formation of a moderately-distorted 
WC at T = 0. 

Results at finite temperature. — At half filling iy = 1/2) 
the ground state of the clean system (W = 0) is a WC with 
a bipartite charge pattern. For a WC the DOS is expected to 
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FIG. 1: (Color online) Top: DOS for the 3D CG model for W = 
0.20. The data are well fit by p{E) ^ E\ S ^ 2 (dashed lines ai'e 
a guide to the eye) around the Fermi level. Center: Same as in the 
top panel for W = 0.40. The insets show the full DOS. Both panels 
have the same horizontal range. Bottom: DOS for the 3D RD model. 
For all W studied the data show a bimodal structure with peaks at 
\E\ ~ 1 and a hard gap of size ~ 2, in stark contrast to the CG 
model. 
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TABLE I: Top: Simulation parameters for tlie simulations of the 
3D Coulomb glass model with Gaussian disorder of strength W at 
finite temperature. L is the system size, TVsa is the number of disor- 
der samples, A^sw is the number of equilibration sweeps, Tmin is the 
lowest temperature, Tmax ~ 0.455 the highest temperature and A^r 
the number temperatures used in the exchange Monte Carlo method. 
Temperatures are measured in units of / Ka, a — 1 being the lattice 
constant. Bottom: Parameters for the 3D RD model simulations. 
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be two delta functions, separated by a charge gap i^wc- The 
energy required to move a particle from a site on the occupied 
sublattice to a site on the unoccupied sublattice is i^wc ~ 2 in 
units of / Ka. Since the peaks of the DOS of the CG are ap- 
proximately centered around \E\ ^ 1 (Fig.[Tl inset), it needs 
to be verified that the observed DOS is indeed representative 
of a glassy phase and not only that of a distorted WC. There- 
fore we study the nature of the phase at finite T by computing 
both an order parameter for a glassy state. 



N 



9gl 



N 



^«-l/2)(nf-l/2), 



i=i 



and an order parameter for the competing Wigner crystal 

2 ^ 

i=l 



(2) 



(3) 



In Eq. ^ a and (3 refer to two copies of the system with 
the same disorder 14511 . If the system forms a Wigner crys- 
tal we expect [(mwc)]av ^ 1 for T < Tc, whereas if the 
system freezes into a glass we expect [((7GL)]av ^ 1 and 
[(mwc)]av for T ^ 0. Here [• • • ]av denotes the average 
over disorder and (• • • ) is a thermal average. 

To locate the putative glass transition we compute the two- 
point finite-size correlation length |460 given by 



1 



2sin(|k„iin|/2) 



x(o) 

X(kmin) 



- 1 



1/2 



(4) 



where kmin = (27r/L, 0, 0) is the smallest nonzero wave vec- 
tor and x(k) is the Fourier transform of the susceptibility 
X = [(^gl) ^ ('?GL)^]av We use four replicas to compute 
X to avoid biases. Because ^/L ^ X[L^/^{T — Tc)], a phase 
transition at Tc is signaled by the correlation lengths for dif- 
ferent L's crossing at the same T = Tc- 

Figure |2] (top panel) shows the qQi^{T) and ni^^{T) or- 
der parameters as a function of temperature for different dis- 
order strengths in the CG model. The glass order parameter 
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no. 2: (Color online) Top: Wigner crystal order parameter 
[(m^c)]av and glass order parameter [(gGL)]av as a function of tem- 
perature T for different disorder strengths W for the CG model. In 
all cases [(mwc)]av <C [(gGL)]av. Center: Finite-size correlation 
length as a function of T for different disorder strengths and sys- 
tem sizes for the CG model. The data show no crossing, i.e., the 
absence of a thermodynamic transition for the studied temperature 
range. Bottom: Wigner crystal order parameter for the RD model 
(L = 8) as a function of temperature for different W. For T < 0.1, 
which Quantitatively agrees with the critical temperatures estimated 
in Ref. I3al , crystalline order emerges. 
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increases as the T 0, whereas the Wigner crystal order pa- 
rameter does not exhibit any ordering tendency, myfc(T) re- 
maining ~ 40 (140) times smaller than qohiT) for W = 0.2 
(W = 0.4) at T = 0.08. Figure |2] (center panel) shows 
the correlation length for the glass order parameter as a func- 
tion of T for the CG model. The data do not cross for the 
studied temperatures and thus there is no sign of a transi- 
tion for T > 0.03, disagreeing with mean-field predictions 



|25t|26t]. The lack of a transition is mirrored by the small cor- 
relation length and the proximity to the ground-state energy 
(not shown). 

In Fig. 12] (bottom panel) we show m^^-, for the RD model 
for disorder strengths up to W = 0.8 covering the disorder 
range studied in Ref. 113811 . For all W studied, rises 
noticeably (in contrast to the CG model). This further under- 
lines that — for the studied disorder range — the phase transi- 
tion in the RD model occurs into a surprisingly robust dis- 
torted Wigner crystal phase. 

Conclusions. — We have analyzed the Coulomb glass at 
low and zero temperature and find that the gap exponent of the 
density of states is close to 5 w 13 — 1 in 3D systems. Further- 
more, we find no evidence of a finite-temperature transition 
into a CG phase in 3D for W — 0.2 and 0.4. This suggests that 
the CG in 3D is at or below its lower critical dimension, which 
would explain the discrepancy with the mean-field results pre- 
dicting a finite transition temperature. Finally, we have shown 
that in a broad disorder range the random-displacement ver- 
sion of the CG model orders into a distorted Wigner crystal 
and not into a glassy state. 
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